#import libraries
library(psych)
library(dplyr)
library(pwr)
library(pscl)
library(car)
library(effects)
library(effectsize)
library(effsize)
library(interactions)
library(jtools)



#load data
data <- read.csv("WP Additional Moderation Study.csv")
View(data)



#recode categorical vars with values
data$gender_coded <- recode_factor(data$gender, 
                                   '1' = 'male', 
                                   '2' = 'female',
                                   '3' = 'non-binary',
                                   '5' = 'prefer not to say',
                                   .default = 'missing')

data$inequality_coded <- recode_factor(data$inequality,
                                       '-1' = 'low',
                                       '1' = 'high')

data$exp_specific_coded <- recode_factor(data$exp_specific,
                                         '1' = 'experienced',
                                         '2' = 'inexperienced')

data$exp_general_coded <- recode_factor(data$exp_general,
                                        '1' = 'experienced',
                                        '2' = 'inexperienced')



#get the sample information
##condition
inequality_count<-table(data$inequality_coded) #low = 81 vs. high = 80
inequality_count

##gender
gender_count<-table(data$gender_coded) #male = 73 vs. female = 76, 6 non-binary, 4 prefer not to say, 2 missing
gender_count

##age
summary(data$age) #41.48
sd(data$age) #12.11

##exp
table(data$exp_specific) #w/exp = 12 vs. w/o exp = 149
table(data$exp_general) #w/exp = 9 vs. w/o exp = 152

#manipulation check
##calculate the manipulation check of inequality 
data <- data %>%
  mutate(manipulation_check = (check_1 + check_2)/2)


##reliability 
inequality_alpha <- data[, c('check_1', 'check_2')]
cronbach.alpha(inequality_alpha) #.90
cor(inequality_alpha) #.83


##check
leveneTest(manipulation_check ~ inequality_coded, data = data) #p = .07
t.test(manipulation_check ~ inequality_coded, data = data, 
       alternative = "two.sided", var.equal = TRUE) #t = -16.35, df = 159, p < .001
cohen.d(manipulation_check ~ inequality_coded,  data = data)

aggregate(manipulation_check ~ inequality_coded, data, mean) #low = 3.01 vs. high = 5.82
aggregate(manipulation_check ~ inequality_coded, data, sd) #low = 1.01 vs. high = 1.16


#regression
##descriptive - willingness
aggregate(willingness ~ inequality_coded, data, mean) #low = 5.57 vs. high = 5.16
aggregate(willingness ~ inequality_coded, data, sd) #low = 1.16 vs. high = 1.27

##descriptive - similarity
aggregate(similarity ~ inequality_coded, data, mean) #low = 4.11 vs. high = 4.03
aggregate(similarity ~ inequality_coded, data, sd) #low = 1.36 vs. high = 1.47


##mean-centering similarity
mean(data$similarity) #Mean = 4.07
sd(data$similarity) #SD = 1.41
median(data$similarity) #Median = 4

data <- data %>%
  mutate(similarity_mc = similarity - mean(similarity))

model = lm(willingness ~ inequality * similarity_mc, data = data) 
summary(model)
standardise(model)

##JN analysis
JN <- johnson_neyman(
  model,
  pred  = inequality,   # focal predictor
  modx  = similarity_mc,   # moderator
  alpha = 0.05 
)
print(JN)

jn$bounds #[-25.01309979,   0.09616524 ]
JN_point <- mean(data$similarity)+0.09616524 
JN_point

percent_below <- mean(data$similarity < JN_point) * 100
percent_below #56.52174
